function [gbest,gbestval,fitcount,fitness]= Template_func(fhd,Dimension,PopulationSize,MaxIter,lowerbound,upperbound,pos,bool_Plot,varargin)

rand('state',sum(100*clock));
mi = MaxIter;
ps = PopulationSize;
D = Dimension;
[pub, plb, vub, vlb] = SetBound(D,ps,lowerbound,upperbound);

%  Parameters of the Algorithm //////////////////////////////////////
cc=[2 2];  % For PSO
iwt=0.9-(1:mi).*(0.5./mi); % For PSO
G0 = 0.28; % For SGSA
epsilon = 1e-3; % For SGSA
g = 0.13; % For QGSA
NumSC = 6; % For Social Class
%  Parameters of the Algorithm //////////////////////////////////////

%  Settings of the Algorithm //////////////////////////////////////
size_each_level = repmat(floor(ps/NumSC),1,NumSC);
size_each_level(1) = ps - (size_each_level(1)*(NumSC-1));
lb_level(1,1) = 1;
ub_level(1,1) = size_each_level(1);
for l = 2 : NumSC
    lb_level(l,1) = sum(size_each_level(1,1:l-1))+1;
    ub_level(l,1) = sum(size_each_level(1,1:l));
end
%  Settings of the Algorithm //////////////////////////////////////

vel = zeros(ps,D); % initialize the velocity of the particles
acc = zeros(ps,D); % initialize the acceleration of the particles

fitness = feval(fhd,pos', varargin{:});
fitcount = ps;
pbest = pos;
pbestval = fitness;
[gbest_fit, gbest_id] = min(pbestval);
gbest = pbest(gbest_id, :);
gbestrep = repmat(gbest, ps, 1);
cbest_fit = min(fitness); % Current best
cworst_fit = max(fitness); % Current worst
gbestval = gbest_fit;

if(bool_Plot)
    h = figure;
    haxes = plot( 0 , 0 );
    XArray = [1];
    YArray = [gbestval];
    title('TEMPLATE');
end

for iter = 2 : mi
    if(cbest_fit == cworst_fit)
        mass = repmat((1/ps),ps,1);
    else
        mass = (fitness'- cworst_fit) / (cbest_fit-cworst_fit);
        mass = mass / sum(mass);
    end
    
    acc = cc(1).*rand(ps,D).*(pbest-pos) + cc(2).*rand(ps,D).*(gbestrep-pos);
    vel = iwt(iter).*vel + acc;
    vel = ((vel>=vlb)&(vel<=vub)).*vel+(vel<vlb).*vlb+(vel>vub).*vub;
    pos = pos+vel;
    pos = ((pos>=plb)&(pos<=pub)).*pos+(pos<plb).*plb+(pos>pub).*pub;
    
    fitness = feval(fhd,pos', varargin{:});
    fitcount = fitcount+ps;
    tmp = (pbestval < fitness);
    temp=repmat(tmp',1,D);
    pbest=temp.*pbest+(1-temp).*pos;
    pbestval=tmp.*pbestval+(1-tmp).*fitness; % update the pbest    
    [gbest_fit, gbest_id] = min(pbestval);
    gbest = pbest(gbest_id, :);
    gbestrep = repmat(gbest, ps, 1);
    cbest_fit = min(fitness); % Current best
    cworst_fit = max(fitness); % Current worst
    gbestval = gbest_fit;
        
    if(bool_Plot)
        XArray = [ XArray iter]; 
        YArray = [ YArray gbestval];
        set( haxes , 'XData' , XArray , 'YData' , YArray );
        drawnow
    end
end

end


